! from "Numerical Recipes"
! hash two integers in two integers. it is strongly assumed that the size
! of the default integer is 4 byte
! used to start seeds for each mpi_rank e for each random index
SUBROUTINE psdes_safe_s(lword,rword)

! USE nrtype

IMPLICIT NONE

INTEGER, INTENT(INOUT) :: lword,rword
INTEGER, PARAMETER :: NITER=4
INTEGER, DIMENSION(4), SAVE :: C1,C2
DATA C1 /Z'BAA96887',Z'1E17D32C',Z'03BCDC3C',Z'0F33D1B2'/
DATA C2 /Z'4B0F3B58',Z'E874F0C3',Z'6955C5A6',Z'55A7CA46'/
INTEGER :: i,ia,ib,iswap
REAL(kind(1.d0)) :: alo,ahi
do i=1,NITER
iswap=rword
ia=ieor(rword,C1(i))
alo=real(iand(ia,65535),kind(1.d0))
ahi=real(iand(ishft(ia,-16),65535),kind(1.d0))
ib=modint(alo*alo+real(not(modint(ahi*ahi)),kind(1.d0)))
ia=ior(ishft(ib,16),iand(ishft(ib,-16),65535))
rword=ieor(lword,modint(real(ieor(C2(i),ia),kind(1.d0))+alo*ahi))
lword=iswap
end do
CONTAINS
!BL
FUNCTION modint(x)
REAL(kind(1.d0)), INTENT(IN) :: x
INTEGER :: modint
REAL(kind(1.d0)) :: a
REAL(kind(1.d0)), PARAMETER :: big=huge(modint), base=big+big+2.d0
a=modulo(x,base)
if (a > big) a=a-base
modint=nint(a,kind=kind(1))
END FUNCTION modint

END SUBROUTINE psdes_safe_s

! program test_psdes
! 
! integer :: il,ir
! 
! il = 2  ;  ir = 3
! print*,'befoer il,ir: ',il,ir
! call psdes_safe_s(il,ir)
! print*,'after il,ir: ',il,ir
! 
! il = 3  ;  ir = 2
! print*,'before il,ir: ',il,ir
! call psdes_safe_s(il,ir)
! print*,'after il,ir: ',il,ir
! 
! end program test_psdes
